library(ggplot2)
library(tidyr)
library(estimatr)
library(tidyverse)
library(lmtest)
library(lsmeans)


# Figure 7: Comparing slopes of eCDFs

setwd("/Users/dmin/Dropbox/Replication_Data/Paper_3/")

Dat <- read.csv("Paper_3_Data_Clean.csv", stringsAsFactors = FALSE)


x =seq(0, 100, by = 1)

my_ecdf <- 
  . %>% 
  split(.$NuclearTheory_fac) %>% 
  map(~ecdf(.$outcome_custom)(x)) %>% 
  data.frame() %>% 
  mutate(x = x) %>% 
  pivot_longer(-x) %>% 
  mutate(
    NuclearTheory_fac = name,
    outcome_custom = value
  )


gg_df <- 
  Dat %>% my_ecdf()


gg_df$name <- as.factor(gg_df$name)

m.interaction <- lm(value~x*name, data=gg_df)
anova(m.interaction)
m.lst <- lstrends(m.interaction, "name", var="x")
m.lst
pairs(m.lst)


g<-ggplot(gg_df,
       aes(y = outcome_custom, x = x, group = NuclearTheory_fac, color = NuclearTheory_fac)) +
  geom_point(alpha=0.4) + geom_smooth(method="lm") + 
  ggtitle("") + xlab("Likelihood of Success") + ylab("Proportion Supporting Use") +
  theme_bw() +
  scale_color_manual(values = c(
    "yellow3",
    "green3",
    "orange3",
    "red3")) +
    annotate(geom = "label", x = 25, y = .13, label = "B = 0.0047", color = "green3") +
    annotate(geom = "label", x = 25, y = .26, label = "B = 0.0075", color = "orange3") +  
    annotate(geom = "label", x = 25, y = .39, label = "B = 0.0066", color = "yellow3") +
    annotate(geom = "label", x = 25, y = .52, label = "B = 0.0062", color = "red3") +
    coord_cartesian(xlim=c(0,100)) + 
    theme(plot.title = element_text(hjust = 0.5), 
          plot.subtitle = element_text(hjust = 0.5),
          legend.position = "none")

g
ggsave("Image_7.png", plot=g, width=2000, height=2000, units=c("px"))
